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Abstract 

We investigate the fragmentation dynamics of an atomic chain under tensile stress. We have 
classified the location, stability type (indices) and energy of all equilibria for the general n-particle 
chain, and have highlighted the importance of saddle points with index > 1. 

We show that for an n = 2-particle chain under tensile stress the index 2 saddle plays a central 
role in organizing the dynamics. We apply normal form theory to analyze phase space structure 
and dynamics in a neighborhood of the index 2 saddle. We define a phase dividing surface (DS) 
that enables us to classify trajectories passing through a neighborhood of the saddle point using 
the values of the integrals associated with the normal form. We also generalize our definition of 
the dividing surface and define an extended dividing surface (EDS), which is used to sample and 
classify all trajectories that pass through a phase space neighborhood of the index 2 saddle at total 
energies less than that of the saddle. 

Classical trajectory simulations are used to study single versus double bond breakage for the 
n = 2 chain under tension. Initial conditions for trajectories are obtained by sampling the EDS at 
constant energy. We sample trajectories at fixed energies both above and below the energy of the 
saddle. The fate of trajectories (single versus double bond breakage) is explored as a function of 
the location of the initial condition on the EDS, and a connection made to the work of Chesnavich 
on collision-induced dissociation. A significant finding is that we can readily identify trajectories 
that exhibit bond healing. Such trajectories pass outside the nominal (index 1) transition state 
for single bond dissociation, but return to the potential well region, possibly several times, before 
ultimately dissociating. 

PACS numbers: 05.45.-a, 31.15.-p, 34.10,+x, 36.20.-r, 62.25.-g 
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I. INTRODUCTION 



There has been much recent interest in the rapidly developing field of 'mechanochem- 
istry', where applied force (e.g., tensile stress) is employed to alter absolute rates and/or 
product ratios of chemical reactions [THE]- A fundamental understanding of the intramolec- 
ular dynamics and reaction kinetics of molecules subject to a tensile force [6-11 J is needed 
to provide a solid theoretical foundation for mechanochemistry, as well as for theories of 
material failure under stress [121411] . polymer rupture [12 ] [TBI [T5T422] adhesion [23], friction 
[21], and biological applications of dynamical force microscopy [2"5T - |3"T] . 

Much previous work on the dynamical consequences of the application of tensile stress 
has focussed on the investigation of fragmentation kinetics of linear chains. Studies of energy 
transfer and equipartition in single chains of coupled anharmonic oscillators have played an 
essential role in the development of nonlinear dynamics, beginning with the seminal work 
of Fermi, Pasta and Ulam [321 [33] (see, for example, refs IM1 - PE6]) . 

The dissociation of a 1-D chain subject to constant tensile force is a problem in unimolec- 
ular kinetics, and a fundamental issue in unimolecular kinetics concerns the applicability of 
statistical approaches such as RRKM [17H5T] or transition state theory [52H57] . Previous 
theoretical work has suggested that dissociation of atomic chains under stress is not amenable 
to simple statistical approaches jlO | I58H66] . Early trajectory simulations on the dynamics of 
Morse chains [HEMD] showed that simple bond stretching or force criteria for bond rupture 
were inadequate, in that apparently broken bonds were observed to reform (bond healing). 
Subsequent simulations of the fragmentation of 1-D Lennard- Jones (LJ) chains at constant 
strain with inclusion of a frictional damping term and a stochastic force modelling inter- 
action with a heat bath showed that healing of incipient breaks is highly efficient [63j (see 
also refs ISTi 68]) . Nonexponential decay, failure of RRKM theory, and extensive transition 
state recrossing effects were found by Bolton, Nordholm and Schranz in their studies of the 
dissociation of 1-D Morse chains (N = 2 — 20) under stress [40J . 

Standard harmonic classical TST has been applied to the dissociation of a 1-D Morse chain 
[HU[65], with the transition state for dissociation of a given bond located at the maximum of 
the effective potential (see below). The harmonic canonical TST rate constant did not agree 
with molecular dynamics calculations, but effects of anharmonicity [101151] on the predictions 
of TST were not systematically investigated. Both RRKM (fully anharmonic, Monte Carlo) 
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and RRK (harmonic appproximation) theory were applied to predict bond dissociation rate 
constants as a function of energy and tensile force for Morse chains under tensile stress [69J . 
For chains with N > 3 atoms a hybrid statistical theory was used involving a harmonic 
approximation for motion at the transition state for bond dissociation [69]. More recent 
work has examined isomerization dynamics for a Morse chain under constant strain subject 
to periodic boundary conditions [TO] . 

Transition state theory, which has long been a cornerstone of the theory of chemical 
reaction rates [5"2TfoT] . has been the subject of renewed interest in recent years [7THS2]- It 
has been established both theoretically and computationally that index one saddles [93] of 
the potential energy surface [HU [95] give rise to a variety of geometrical structures in phase 
space, enabling the realization of Wigner's vision of a transition state theory constructed in 
phase space [7H492] . 

Following these studies, attention has naturally focussed on phase space structures as- 
sociated with saddles of index greater than one, and their possible dynamical significance 
[91)H99] . We have investigated phase space structures and their influence on transport in 
phase space associated with index two saddles of the potential energy surface for n degree- 
of-freedom (DoF) deterministic, time-independent Hamiltonian systems [96J [99] . We have 
shown that, for isomerization dynamics in a model n = 2 DoF potential, it is possible to 
distinguish between 'concerted' and 'stepwise' (or sequential) isomerization trajectories in 
a dynamically significant way using a normal form Hamiltonian describing the phase space 
structure in the vicinity of the index 2 saddle. The importance of index 2 saddles for the 
question of dynamical separability of tight versus roaming mechanisms has recently been 
established [lOOj . As discussed in detail in the present paper, higher index saddles arise 
naturally in the problem of an atomic chain under tensile stress, and we are able to apply 
the methods and insights from our previous studies to this problem. 

In the present work normal form theory together with classical trajectory simulations are 
used to investigate the fragmentation kinetics and phase space structure of short tethered 
atomic chains under constant tensile stress. Following previous work, we model the inter- 
atomic interactions using Morse potentials [69]. Most of our work concerns a 'chain' with 
n = 2 particles. Our focus is on the relation between the phase space structure in the vicin- 
ity of the index 2 saddle in this system, as described by a normal form Hamiltonian, and 
the fate of trajectories passing through a suitably defined dividing surface. The competition 
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between single and double bond breaking is of particular interest, as is the phenomenon of 
bond 'healing' Ml EE] ■ 

Our phase space approach to the problem of single versus double bond breakage in the 
Morse chain under tensile stress is related to the work of Chesnavich on collision induced 
dissociation |1011 1102] ; the results obtained in the present paper serve to confirm and extend 
Chesnavich's insights concerning the nature of the boundaries between phase space regions 
associated with single and double bond breakage. 

The structure of this paper is as follows. In Sec. [TT] we describe the Hamiltonian for a 
tethered Morse oscillator chain under tension. We analyze the location and types (indices) 
of equilibria for the general n-particle chain, and provide numerical results for the n = 2 
case. Specifically, we note the presence of an index 2 saddle in the 2-atom chain under tensile 



stress. In Section [HI] we discuss phase space structure and dynamics in a neighborhood of a 
saddle-saddle equilibrium point. Based on previous work [96j EH], we define a phase space 
dividing surface (DS) that enables us to classify trajectories passing through a neighborhood 
of the saddle point using the values of the integrals associated with the normal form. We 
also generalize our definition of the dividing surface and define an extended dividing surface 
(EDS); the EDS can be used to sample and classify trajectories that pass through a phase 
space neighborhood of the index 2 saddle at total energies less or greater than that of the 



saddle. In Section[[V]we describe our classical trajectory studies of single versus double bond 
breakage for the n = 2 chain under tension. Initial conditions for trajectories are obtained 
by sampling the EDS at constant energy, and we sample trajectories at fixed energies both 
above and below the energy of the saddle. The fate of trajectories (single versus double bond 
breakage) is explored as a function of the location of the initial condition on the EDS. We 
find that, at the boundary between regions of the EDS associated with trajectories exhibiting 
breakage of one or the other bond, there is either (i) a region coresponding to double bond 
breakage, or (ii) a trajectory that is trapped in the vicinity of the potential minimum. 
The existence of a layer of doubly dissociative trajectories corresponds to the 'coating' 
phenomenon noted by Chesnavich in his work on collision-induced dissociation |101l 1102] . 
Our trajectory results also enable us to identify trajectories that exhibit bond healing; such 
trajectories pass beyond the nominal (index 1) transition state for single bond dissociation, 
but return to the well region possibly several times before ultimately dissociating. Section [V] 
concludes. The procedure for sampling the DS and EDS used in the computations reported 
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here is described in detail in Appendix [Aj 
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II. MORSE OSCILLATOR CHAIN UNDER TENSION: HAMILTONIAN AND 
EQUILIBRIA 

In this section we describe the Hamiltonian used to model a ID chain of particles under 
tensile stress. After deriving the equations of motion, we determine the location of equilib- 
rium points for chains composed of different numbers of particles. We first review the single 
particle 'chain' and then give the general form of the Hamiltonian for the n-particle case 
together with analysis of equilibria. 



A. One particle Morse chain 

We consider first a one particle 'chain'. Although this case is trivial, the results obtained 
will be useful when dealing with the general case of n-particle chain. 

The single degree of freedom (DoF) Hamiltonian for our 1-particle chain is the sum of a 
kinetic energy term and a potential energy described by a Morse oscillator potential [103J 
representing a particle of mass m tethered to a wall of infinite mass. To the Morse potential 
we add a linear term in the bond coordinate representing a tensile force exerted on the 
particle, whose magnitude is determined by a parameter /. The one DoF Hamiltonian takes 
the form: 

H(x,p; /) = !- + V M (x) - f(x - x e ) (2.1) 
2m 

where the Morse potential Vm(x) is 

Vm(x) = D [l - exp(-/3(x - x e ))f. (2.2) 

In the following we will measure energies in units of -Do, the unperturbed Morse dissoci- 
ation energy, and length in units of x e , the equilibrium bond distance for zero tensile force. 
The Morse parameter (3 will be set equal to l/x e so that in our units its numerical value is 
P = l. 

Associated Hamiltonian equations of motion (vector field) are: 

* = = ^ (2.3a) 
op m 

OH 

P = = 2D /3{exp [-2p(x - x e )} - exp [-p{x - x e )}} + f. (2.3b) 
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The equilibrium points of this vector field satisfy the following equations 

x{x,p) = (2.4a) 
p(x,p) = 0. (2.4b) 



Eq. (2.4a) is always satisfied for p = 0. The number of (real) roots of eq. (2.4b) depends on 
the value of the / parameter. The critical value for this parameter is f cr n = D (3/2, which 
has the value f crit = \ in our units. For < / < f crit we have two real roots, for / = f crit 
we have one real root and for / > f crit there is no real root. For < / < f crit , the two roots 
are: 

'D Q p+y/DlP-2D Q pf 




2D (3 



+ x e (2.5a) 



, = .: H ^^M +Il (2 . 5b ) 

The equilibrium point (x+,0) is the equilibrium phase space point for the Morse oscillator 
under tensile stress while the (x_,0) solution is a saddle point. (The stability of these two 
solutions can easily be verified by computing the eigenvalues of the matrix associated with 
the linearization of Hamilton's equations about the equilibrium point of interest.) 

For / parameter equal to f crit the two equilibrium points (x + ,0) and (x_,0) merge at a 
single point, which we will denote by (x*,0), in a saddle node bifurcation. For this critical 
value of / the only real root is 

x* = -^ln(2)+x e (2.6) 

As the force parameter / — > 0, x + — > x e , the unperturbed Morse equilibrium distance, while 
,x_ — >■ +oo. Figure [l](a) shows the shape of the 1-D potential for values of / parameter 
between and / cr jt and Figure [T|(b) shows the locations of the different equilibrium points 
discussed above as a function of /. 



B. n-particle Morse chain under tension 

We now consider a chain of n particles tethered to an infinite mass wall. For simplicity 
we assume all particles have the same mass m. 

We use two different sets of phase space coordinates to describe the dynamics. The first 
set of coordinates are 'external' or lab-fixed' variables; the configuration space coordinates 
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are denoted (xi,xz, . . . ,x n ) with conjugate momenta (p xl ,Px 2 , ■ ■ ■ ,Px„)- The second set of 
coordinates are 'internal' or bond coordinates; coordinates are denoted r = (ri,r 2 ,. . . ,r n ) 
with conjugate momenta p r = (p ri ,Pr 2 , ■ ■ ■ ,Pr„)- 

Figure [2] shows the definitions of these two coordinate systems for the case of a 2-particle 
chain. The general relation between coordinates (x\, X2, ■ ■ ■ , x n ) and (7*1, r 2 , . . . , r n ) is 

k 

x k = ^Jrj, k = l,...,n. (2.7) 

3=1 

The kinetic energy in external coordinates is a diagonal sum of quadratic terms. The 
potential energy consists of a sum of pairwise Morse interactions between adjacent particles 
together with potential term linear in the coordinate x n of the last particle (equivalent to a 
constant force applied to this particle). The Hamiltonian is therefore 

n 2 



H(x 1 ,...,x n ,p xl ,...,p Xn ;f) = ^2^- + V(x 1 ,...,x n ;f) 

i=i 

with potential term 

n 

V(x 1 , . . . , x n ; /) = ^ V M (r-i) - f(x n - nx e ) 
i=i 

n 
i=l 

In terms of bond coordinates, the Hamiltonian is (cf. ref. | 



(2.8) 



(2.9a) 
(2.9b) 



H(r,p r ;f) 



n 2 n—1 

2m z — ' m L — ' m 

i=2 k=l 



Pn 



(2.10) 



The vector field associated with Hamiltonian (2.10) is: 



OH 


_ Pn _ 


_ Pr2 


dp ri 


m 


m 


OH 


__ 2Pr t 


Pn-i 


9p n 


m 


m 


dH 


_ 2Pr n 




9pr n 


m 


m 



m 



, n — 1 



Pr 



dH 



2D /3 [exp{-2/3(r 4 - x e )} - exp{-[3( ri - x e )}\ +f,i = l, 



(2.11a) 
(2.11b) 
(2.11c) 
,n. (2. lid) 



The equilibrium points are found by setting the time derivatives in eq. (2.11) to zero. The 



solutions are in fact easily found in terms of equilibrium points for the n—1 DoF system. 
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TABLE I: Pascal's triangle structure of the indices of the different saddles (k) depending on the 
number of particles in the chain (n). Index means a stable centre (well). 



The time derivatives (2.11a), (2.11b) and (2.11c) are zero only for p r . = 0. The solutions 



obtained by setting the time derivatives of the momenta to zero depend on the value of / 
parameter. For < / < / cr jt we have two real roots for each which are related with the 
one DoF roots by r{± = x±, so that the number of equilibria for an n-particle chain is 2 n . 
The stability of these equilibria is determined by the stability of the one DoF equilibrium 
points, so that we obtain a c-centre-s-saddle when the x + root occurs c times and the X- 
root s times with n = c+s. The organisation of the different saddle indices follows a Pascal's 
triangle structure as shown in table [T| 

In the rest of this paper we focus on the n = 2-particle chain. For this case we have four 
equilibrium points for < / < fcrit and only one for / = f crit . The situation is summarized 
as follows: 

• < / < f crit : EPx = (r 1+ ,r 2+ ,0,0), EP 2 = (r 1+ , r 2 _, 0, 0), EP 3 = (r^, r 2+ , 0, 0), 
EP 4 = (n_,r2_,0,0). 

• / = fcrit- EP * = (ri*,r 2!M 0,0). 

EPi is a center, EP 2)3 are index 1 saddles and EP 4 is an index 2 saddle. Figure [31(a) shows 
contours of the potential function in (r 1; r 2 ) space as well as the location of the equilibria 
while Figure (3^b) shows the evolution of the energy of these equilibria points as a function 
of the force parameter /. 
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III. PHASE SPACE STRUCTURE AND DYNAMICS IN A NEIGHBORHOOD 
OF THE INDEX 2 SADDLE 

In the present work we focus on the role of the index-2 saddle EP4 in the bond dissociation 
kinetics of the 2-particle chain. We study fragmentation dynamics at energies close to the 
energy of this equilibrium point. Our aim here is to investigate the competition between 
single and double bond breakage for trajectories that enter the phase space neighborhood 
of the index 2 saddle. We will study trajectories initiated on a neighbourhood of EP 4 and 
investigate the dynamical role of the index 2 saddle on the fragmentation of the chain. 

In order to carry out this program we construct a normal form Hamiltonian which pro- 
vides an integrable approximation to the dynamics associated with the original Hamiltonian 
in a neighbourhood of the saddle point EP4. In this section we first briefly describe the 
construction of a normal form Hamiltonian in a neighbourhood of a saddle point; we then 
show how this Hamiltonian can be used to describe the dynamics in this neighbourhood 
and different phase space objects such as normally hyperbolic invariant manifolds (NHIMs), 
dividing surfaces and extended dividing surfaces (EDS). 

A. The normal form 

Normal form theory has been mainly applied in reaction dynamics involving equilibrium 
points of saddle x center x • ■ ■ x center stability type (that is, index 1 saddles). For this kind 
of saddle the normal form Hamiltonian permits the construction of certain phase space 
structures which are of central importance for the reaction dynamics (for a review and 
further references, see ref. 80). 

However normal form theory is also useful for describing dynamics in the vicinity of equi- 
libria of different stability type such as EP4 in the present problem, which is of saddle x saddle 
stability type. Recent work has begun to investigate reaction dynamics mediated by higher 
index saddles [96H99]. 

The construction of the normal form coordinate set is carried out by Poincare-Birkhoff 
normal form theory. The details of this theory are by now well-known and have been 
presented in detail in a number of reviews and books (see, for example, refs l80l 11041) . Normal 
form theory provides an algorithmic procedure for finding a non-linear symplectic change of 
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variables, 



(q,p) = T(x,p x ), 



(3.1) 



which will transform a given Hamiltonian into a new, simpler, Hamiltonian, 



H NF (q,p) = HiT- 1 (q,p)) = tf(x,p x ). 



(3.2) 



The form of the resulting Hamiltonian H^f is constrained by imposing conditions such 
as the requirement that H^p should Poisson commute with a certain Hamiltonian H Q ; in 
such a case we say that the resulting Hamiltonian is in normal form with respect to this 
Hamiltonian Hq. The Hamiltonian Hq is often taken to have the simplest possible form, 
namely, a Hamiltonian describing n uncoupled harmonic oscillators. 

The normal form Hamiltonian will describe in a neighbourhood C of the equilibrium point 
an integrable system which decouples the dynamics into "reaction coordinates" and "bath 
modes" . 

In the general case of a n DoF Hamiltonian system, the matrix associated with the 
linearization of Hamilton's equations about the equilibrium point has k pairs of real eigen- 
values of equal magnitude, but opposite in sign (±Aj) and n — k pairs of complex con- 
jugate purely imaginary eigenvalues (±Wj). We will assume that a non-resonance condi- 
tion holds between the eigenvalues, i.e., the eigenvalues (uJk+i, ■ ■ ■ ,ou n ) satisfy the relation 
Ck+iUJk+i + • • ■ + c n u n 7^ for any vector of integers (ck+i, ■ ■ ■ , c n ). Under these conditions the 
normalisation procedure transforms the original Hamiltonian to an even order polynomial 
in the normal form coordinates (qi, . . . , q n ,Pi, ■ ■ ■ ,Pn)- 

The normal form Hamiltonian thus obtained, H NF , describes an integrable system which 
approximates the dynamics of the original Hamiltonian H in a neighbourhood of the equi- 
librium point. As the resulting Hamiltonian is integrable, we can find n integrals of motion 
and express the normal form Hamiltonian explicitly in terms of these integrals: 



where the higher order terms (hot) are at least of order two in the integrals (ii, . . . , I n ). The 
expressions of these actions in terms of the normal form coordinates (q±, . . . , q n ,Pi, ■ ■ ■ ,Pn) 



HNF(q,p) = K(I U ..., 




(3.3a) 
(3.3b) 



= Xih H h \ k h + ^k+ih+i H h Wn4 + HOT, 
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are: 



h = qiPi, i = l,...,k (3.4a) 

Ii = \{<i+P)), J = k + l,...,n. (3.4b) 

The vector field associated with the normal form Hamiltonian H^f is 

dH NF dK dli . , 

qi = -v, — = jn~jr = A *qi ( 3 - 5a ) 

dpi dli dpi 

^ = dq^~ = ~~dldq- = ~ AiPu i = 1 >---> k ( 3 - 5b ) 

q.-™™-?IL?L-n. p . (3 5c) 

q '~ dp, ~ dI j dp j -"' P ' [ j 

p j = -^ = -™¥j- = -n jqj , j = k + l,...,n (3.5d) 

where we have defined the frequencies 

dK . 

= ~d~T' l ^ " ^ ^ 

toj = -gf, j = k + l,...,n. (3.6b) 

We now introduce a new set of canonical coordinates, (Qj,Pj), for the saddle planes. These 
coordinates are useful for describing the dynamics in each saddle plane and are obtained by 
a canonical transformation of the NF variables (qi,Pi) i — 1, . . . , k 

Qi = ^=(Qi-Pi) ( 3 - 7a ) 

Pi = -j=(qi+Pi), i = (3.7b) 

In terms of these coordinates the action variables are 

U = \(F?-<%\ < = l,...,fc (3.8) 

while the vector field for the saddle modes transforms to: 

■ dH NF dK dh 

Qi = ^p- = dL t m = AtPi (3 - 9a) 

• dH NF dK dh 

= dQ^ = ~~d~Ti~dQi = ^ h « = 1.-.* (3-9b) 

and the vector field for the bath modes is unchanged. 
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B. Phase space structures in normal form coordinates 



1. Crossing and non-crossing trajectories 

The normalization algorithm provides us with an integrable system. The dynamics sepa- 
rates into independent motions in the position-momentum planes for the saddle modes and 
the bath modes, so that the full dynamics is the cartesian product of the motion in these 
planes. 

There are k saddle ('reactive') planes spanned by coordinates (qi,Pi) or (Q^Pj), i = 
1, . . . , k. The phase space for each saddle mode is foliated by manifolds specified by the 
value of the integral of motion ij, i — 1, . . . , k. For each value of Ii, the trajectory curves 
(solutions of Hamilton's equations) are simply the two branches of the hyperbola, 



The dynamics in the remaining n — k bath modes consists of uniform rotation in angle 6j 
conjugate to the conserved action Ij in the respective (qj,Pj) planes with j — k + 1, . . . , n. 
The total energy E of the normalized system is a function of the action variables only: 



The k saddle integrals Ii, i = 1, . . . , k, can either be positive or negative, leading to two 
different types of branches of the hyperbola in each saddle plane depending on the sign of the 
integral. For Jj > the two branches of the hyperbola form what we call crossing trajectories 
whereas those branches for which Ij < are called non-crossing trajectories. Considering 
constant action hyperbolae in the (Qi, Pi) plane, crossing trajectories are those for which the 
sign of Qi changes along the branch whereas non-crossing trajectories are those for which 
the sign of Qi remains the same. For crossing trajectories, we can distinguish trajectories 
for which the sign of Qi changes from negative to positive from trajectories for which the 
sign of Qi changes from positive to negative. 

As discussed in ref. [HQ we can introduce a symbolic description of trajectories in the 
neighbourhood of the origin in the saddle plane. A trajectory is labelled by 2 symbols, (/; i), 
where i — ±, / = ±. Here i refers to the initial sign of Qj as it enters the neighbourhood of 
the origin and / refers to the final sign of Qj as it leaves the neighbourhood. In the case of 
multiple saddle mode planes, we can extend this symbolic description to 2k indices by taking 




(3.10) 




(3.11) 
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the cartesian product of the k saddle planes and labelling a trajectory by (fi . . . fk', H ■ ■ ■ ik)- 

For example, with two saddle modes, a trajectory labelled (H — ; ) crosses in the first plane 

and does not cross in the second plane, the sign of Q2 remaining negative. Figure [4] shows 
a schematic representation of the different hyperbola in the two saddle planes. 

2. Dividing surfaces, extended dividing surfaces and NHIMs 

The notion of a dividing surface originates in the study of chemical reaction dynamics. In 
this context one is interested in defining a surface in phase space separating reactants from 
products, through which all reactive trajectories must pass, and which is never encountered 
by any nonreactive trajectories. This surface is conventionally referred to as the transition 
state (TS). Transition state theory (TST) has a long history going back to Eyring [105J, 
Wigner |52j and Keck |53j (variational transition state theory, introduced in order to deal 
with the problem of recrossing [55, 56]). 

Many chemical applications of TST employ a TS defined in configuration space [57J. For 
the case of 2 DoF, a dynamically based approach to TST was pioneered by Pechukas, Pollak 
and Child, who introduced the notion of periodic orbit dividing surface (PODS) |54j . 

There has been significant recent progress in generalizing the dynamically based PODS 
approach to obtain a definition of dividing surfaces in phase space for n > 3 DoF. The 
phase space approach uses normal form theory to construct a normal form Hamiltonian 
which reproduces the dynamics in the neighbourhood of the phase space region of interest 
|80j . This normal form Hamiltonian is the key object which enables the precise mathematical 
realization of the intuitive idea of a 'surface of no return' intersecting all reactive trajectories 
jHH]. 

We now turn to the definition of the dividing surface (DS) and the so-called extended 
dividing surface (EDS) for the present problem in which there is a dynamically significant 
index 2 saddle point. The DS and the EDS are codimension one surfaces within the constant 
energy surface that are transverse to the vector field defined by Hamilton's equations, and 
consist of the set of phase space points for which a suitable distance (to be defined) from 
the equilibrium point is (locally) minimal along the trajectory passing through the phase 
space point. 

Whereas the energy at which the DS is defined cannot be arbitrary, the EDS is defined 
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for any value of the energy. For simplicity we take the zero of energy to be the energy of the 
equilibrium point in the neighbourhood of which we compute the normal form Hamiltonian. 
We define the square of the distance from the origin in phase space to be 

k 

D = \T,(Q* + P ?) ( 3 - 12 ) 

i=l 

where the sum is taken over the saddle DoF. (Note that a different definition of the distance 
D was used in ref. HH1 For further discussion of this point, see Appendix [Aj) Requiring that 
points on the DS and the EDS are those for which this distance is minimized implies that 
the time derivative of this distance along the trajectory passing through the phase point 
should vanish: 

k 

D = J2 A iQi P i = °- ( 3 - 13 ) 

i=l 

This relation enables us to find those phase space points within a given energy surface that 
satisfy the (local) minimum distance requirement. 

The EDS is defined by the intersection of the following 2n — 1 dimensional surfaces in 
the n-dimensional phase space: 



Si(Qi,Pi, ■ ■ - ,Qk,Pk,qk+i,Pk+i, - ■ ■ ,q n ,Pn) = K(I u ...,l n ) -E = 0, 

k 

S 2 {Qi,Pi,...,Qk,Pk,qk+i,Pk+i,---,q n ,Pn) = ^AiQiPi = (3.14) 



i=i 



The DS is obtained from (3.14) by restricting to points that lie on crossing trajectories. This 
condition implies that the actions associated with the saddle degrees-of-freedom should be 
positive. On the other hand, points on the EDS are not required to lie on crossing trajectories 
so that all phase space points within a given energy surface satisfying the minimum distance 
requirement belong to the EDS. 

In the simplest and most familar case, a chemical transformation (reaction) involves an 
equilibrium point of saddle x centre x • • • x centre stability type. The constant energy dividing 
surface intersects all reactive trajectories and consists of points for which the phase space 
distance from the saddle point is a minimum. Clearly, for a trajectory to react, the action 
associated with the saddle mode must be positive in order for the trajectory to overcome 
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the barrier, whereas the actions relative to the bath are always positive or zero; the energy 
at which the dividing surface is defined can only be positive. The set of points in a given 
energy surface which minimize the distance from the saddle equilibrium point are then 
located precisely on the line Q± = in the saddle mode plane and in a disk in the bath 
mode planes. When all the energy is in the saddle direction, the saddle action is a maximum 
and the actions relative to the bath modes are all zero so that K(Ii max , 0, . . . , 0) = E > 0. 
Conversely, when all the energy is distributed amongst the bath modes, the action of the 
saddle is zero; such phase points do not belong to the DS and in fact constitute a normally 
hyperbolic invariant manifold, or NHIM [104] . For the case of a 2-DoF system with one 
saddle mode and one bath mode this NHIM is a periodic orbit, the so-called periodic orbit 
dividing surface (PODS) [M]- This PODS formally divides the DS into two disjoint pieces: 
one part for which the reactive trajectories cross the DS from reactants to products (forward 
reaction) and one piece for which they cross the DS from products to reactants (backward 
reaction). Again, for a 2-DoF system, each dividing surface (forward/backwards) is a disk, 
so that the union of the DS and the NHIM is topologically equivalent to a sphere, as 
represented in Figure [5^a), where the NHIM (PODS) is the equator of this sphere. The 
NHIM, represented in green, separates the sphere into two hemispheres, the forward (red) 
and backward (blue) hemispheres. 

The phase space definition of the DS is essential for a fundamental understanding of 
reaction dynamics [Ml [HO] . On the other hand, if we are interested in characterising all 
possible types of dynamics in the vicinity of a saddle equilibrium, the DS as defined above 
does not provide complete information, as (by definition) it intersects only those trajectories 
which react. In order to define an object which will capture the totality of the dynamics 
in the phase space neighborhood of the equilibrium point we must extend the notion of 
a dividing surface by relaxing the restriction to crossing trajectories and so include non- 
crossing trajectories as well. 

This new object is called the extended dividing surface (EDS). Formally the definition of 
this surface is a surface within an energy surface consisting of phase points which (locally) 
minimize the distance from the saddle point along their trajectory, without any restriction 
to crossing trajectories. As a consequence, the action associated with the saddle mode can 
now be either positive or negative. While the DS was only defined for positive energies, the 
EDS is defined for energies either positive or negative. The topology of the EDS depends, 
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however, on whether the energy is positive or negative. 

Consider the case of an index- 1 saddle. For positive values of the saddle action I\ the 
EDS is simply the usual DS [HH]; that is to say, it is a sphere with forward and backward 
hemispheres separated by a NHIM. If we allow the saddle action to be negative, for E > 
the equation K(I\ < 0, 1 2 , . . . , I n ) — E > must be satisfied. In order to compensate for the 
loss of energy in the saddle direction, the actions of the bath modes have to increase. Topo- 
logically the EDS is equivalent to a hyperboloid as represented schematically in Figure [51(b). 




Notice that in this case the EDS is not compact as the action of the saddle mode can become 
arbitrarily large and negative, with a corresponding increase of the bath mode actions to 
conserve energy. Of course this picture is valid only in the neighbourhood of validity of the 
normal form Hamiltonian and will in general break down outside this neighbourhood. It is 
important to realise that the EDS only has meaning within this (bounded) neighbourhood 
of validity for the normal form. 

For negative energies the DS is not a subset of the EDS and this EDS is topologically 
equivalent to a two sheeted hyperboloid represented in Figure [5] (c). 

For an index 1 saddle in an n DoF system, the EDS is of subsidiary importance to the 
DS for understanding phase space structure and reaction dynamics. When we consider the 
case of higher index saddles, however, the EDS becomes an object of central importance. 
For this more complex situation the notion of dividing surface itself becomes problematic 
due to the fact that concepts like 'reaction pathway' are not necessarily well-defined (see 
ref. I1UUI and refs therein). Even if we retain the strict mathematical definition of the DS 
as the surface which intersects all reactive trajectories, the resulting surface gives a rather 
restricted view of the dynamics in the vicinity of the equilibrium point, whereas the EDS is 
an object encoding all information about the dynamics in the vicinity of the equilibrium. 

We now consider in detail the case of an n DoF system with an index 2 saddle. In this 
case we have two saddle planes and n — 2 bath mode planes. For both positive or negative 
energies we can always find a combination of the different actions which satisfy the equation 
K(Ii,l2, . . . ,I n ) = E with Ix t 2 positive or negative and the bath actions always positive. 
Points on the EDS projected onto the saddle plane (Qi, Pi) can therefore be located on any 
quadrant of the saddle plane. 

Fixing a point in one of the saddle planes, (Qi,Pi) say, corresponding phase points on 
the EDS must minimize the distance from the equilibrium point. In Appendix|A]we describe 
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our procedure for locating phase points on the EDS at fixed energy. Specifically, we show 
how a parametrization of the hyperbola describing the NF dynamics in each saddle plane 
can be used to sample the EDS. Parameters t iy i = 1, . . . , k, describe the position of the 
phase point on these hyperbolae; a useful representation of the EDS therefore utilizes the 
space of these ti parameters. 

As discussed in Appendix [A] such a representation is multivalued: for every choice of the 
ti parameters, there will be several phase space points on the EDS corresponding to different 
combinations of signs of the actions I{. It is useful to have a representation of the EDS that 
retains crucial information such as location of a phase point in the saddle planes and the 
distance from the equilibrium. In a saddle plane, the character of a trajectory depends on 
which quadrant the trajectory intersects the EDS. The location of the phase point in the 
saddle plane can be specified by an angle; for the case of an index 2 saddle there are two 
angles. We use these two angles and the distance in phase space from the saddle point to 
construct a toroidal representation of the EDS. Figure [6] shows the definitions of the angles 
used and the toroidal representation of the EDS. The torus topology is appropriate because 
phase space points are periodic functions of the angles ti with period 2n . 

As phase points composing the EDS can belong to crossing or non-crossing trajectories, 
the EDS splits into several parts (four parts for each saddle plane, 4 fc parts for an index 
k saddle). In order to label these different parts of the EDS, we can use the symbolic 
description defined for crossing and non-crossing trajectories. For the case of k = 2 saddle 
modes the 16 parts of the EDS are labelled by four indices, (fi fi\i\ h)- For example, 
the part of the EDS having trajectories which cross in the first saddle plane and which do 
not cross in the second saddle plane (with Q 2 > 0) is denoted EDS( ++; _ + ). Notice that, 
according to our definitions, the DS is a subset of the EDS consisting of the parts with 
symbol codes (++; ), (H — ; — h), ( — h; H — ) and ( ; ++). 

An important question concerns the method used to sample points on the EDS. One 
method was described in a previous work [99]. In the present paper we use another sampling 
method for the EDS, which is described in detail in Appendix [A} 
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IV. SINGLE VERSUS DOUBLE BOND BREAKAGE FOR THE n = 2 CHAIN 



In this section we use classical trajectories to investigate single versus double bond break- 
age in the n = 2 particle chain. By 'double bond breakage', we mean fission of both bonds 
leading to the complete breakup of the chain. Specifically, we investigate the behavior of 
trajectories that pass through a neighborhood of the index 2 saddle. Similar dynamics was 
explored in the pioneering work of Chesnavich on collision induced dissociation reactions 
[TOTl [102] . 

Our problem has two essential parameters. The first is the magnitude / of the tensile 
force. Varying this parameter changes both the absolute and relative energies of the various 
equlibrium points: the well (EPi), the two index one saddles (EP 2 ,3), and the index two 
saddle (EP4). In the present work we set / = 0.1, which represents a physically realistic 
value of the tensile stress |69j. The other parameter of interest is the energy at which we 
study the breakage of the chain. Specifically, we are interested in the energy dependence of 
single versus double breakage of the chain at energies close to the energy of the index two 
saddle (i.e., near threshold). 

To study these questions, we propagate classical trajectories and examine their fate. 
Points on the EDS (defined as in the previous Section) are used as initial conditions for 
trajectory propagation. We examine the dependence of the fate of trajectories on their 
location on the EDS, and attempt to understand the distribution over the EDS in terms 
of single breakage of either the first or the second bond and double breakage. We also 
investigate the evolution of this distribution as we change the energy at which the EDS is 
sampled from an energy greater than that of the index two saddle EP4 to an energy below. 
Finally we discuss a possible manifestation and interpretation of the phenomenon of bond 
healing [63j |67J [68] in the light of the behavior found in our classical trajectory simulations. 

A. General considerations 

To investigate the fate of trajectories at a fixed energy in the neighbourhood of the index 
2 saddle EP 4 for the n = 2 particle chain, we apply normal form theory to construct a 
Hamiltonian which approximates the dynamics in a neighbourhood of the equlibrium EP 4 . 
This equilibrium point is of saddle x saddle stability type (for n = 2 the normal form is a 2 
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DoF Hamiltonian without bath modes). With this NF Hamiltonian in hand, the sampling 
of the EDS is carried out using the procedure described in Appendix |Aj Using the backward 
transformation from normal form coordinates to physical coordinates defined in equations 



(3.1 ) and (3.2), we obtain phase space points belonging to the EDS in the original physical 
coordinates. We use these points as initial conditions for trajectory propagation, and monitor 
the fate of each trajectory. 

Several different types of trajectory are possible. The first, denoted type 0, are trajectories 
which remain trapped in the well region in the neighbourhood of equilibrium (minimum) EPi 
over the full integration time. For these trajectories no bonds rupture. The second type of 
trajectories, denoted type 1, are those which exhibit single bond breakage; those trajectories 
for which bond r\ breaks will be denoted type li and those for which bond r 2 breaks type I2. 
Finally, type 2 trajectories are those for which both bonds break, corresponding to complete 
fragmentation of the chain. 

One important question concerns the criteria used to decide whether or not a bond has 
actually broken. In the work reported here, we propagate trajectories for a long time and 
look at the values of the bond lengths r\ and r 2 at the end of the run. By choosing thresholds 
for both r\ and r 2 coordinates we can determine the trajectory type by comparing the final 
values of the bond coordinates of the propagated trajectory with suitably chosen values. 
For our simulations, we integrate trajectories for 2000 time units (time units as defined in 
Section [il]) and the thresholds in ri and r 2 coordinates were chosen to be 30 (units of r e ). 

The matrix associated with the linearisation of Hamiltons equations about the equilibrium 
point EP4 has two pairs of real eigenvalues of equal magnitude and opposite signs. For the 
present problem, with tensile stress parameter / = 0.1, these eigenvalues are ±Ai = ±0.4972 
and ±A2 = ±0.1899. As is shown in Figure [7], the projections of the associated eigenvectors 
into the (r 1; r 2 ) plane indicate reactive directions in the two saddle planes given by the normal 
form Hamiltonian. Roughly speaking, the 'direction' associated with the first eigenvalue 
provides in a sense an indicator for which bond will be broken. If we cross the saddle in 
this eigendirection from negative to positive Q, we will follow the direction given by the 
eigenvector associated with ±Ai and consequently we are likely to break bond r 2 . If we 
cross in the opposite direction, we follow the direction given by the eigenvector associated 
with — Ai and are likely to break bond r±. The directions corresponding to the eigenvector 
associated with the second eigenvalue indicate whether the trajectory will enter the well 
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region (in the neighbourhood of EPi) or pass into the region where we expect double bond 
breakage to occur. This crude approach to the dynamics for trajectories entering the vicinity 
of the index 2 saddle does not capture the full complexity of the dynamics; in fact there is a 
competition between these two directions in determining the behaviour of the trajectories. 

B. Results 

1. Energy above the index 2 saddle 

We consider first the behaviour of trajectories initiated on the EDS for energies above 
the energy of the EP 4 saddle point. 

To understand our results, we examine the dynamics associated with the second saddle 
plane in normal form coordinates; roughly speaking, motion along the associated eigendi- 
rection determines whether or not the trajectory enters the well region or exits directly into 
a region of configuration space where we might expect double bond breakage to occur. The 
normal form Hamiltonian in this plane has 4 distinct types of dynamics (combinations of 
symbols ±): two correspond to crossing trajectories and two to non-crossing trajectories. 
Non-crossing trajectories remain either in the well region (Q 2 < 0) or outside the well re- 
gion (Q 2 > 0). There are therefore two relevant classes of trajectories. The first class is 
composed of trajectories for which the sign of Q 2 changes from negative to positive together 
with non-crossing trajectories for which Q 2 > 0. The second class is composed of trajecto- 
ries for which the sign of Q 2 changes from positive to negative together with non-crossing 
trajectories for which Q 2 < remain negative. Using the symbolic classification, the first 
class consists of parts of the EDS associated with codes (/i+;«i— ) and (f 1 +;i 1 +), and the 
second class with codes — ; and — ; H-), where f\ and i\ are used to designate all 
the possible symbols for the first saddle plane (Qi, Pi). 

Let us focus on the first of these classes, ) U (/i+;«i+). This class splits into 

two different subclasses: ) and (/i+;ii+). Each of these two subclasses splits in 

turn into several distinct subsubclasses according to the symbols associated with motion in 

the first saddle plane. The first subclass is composed of the parts (++; H — ), (++; ), 

( — h; ) and ( — h; H — ), whereas the second subclass is composed of the parts (++; ++), 

(++; — h), ( — h; — h) and ( — h; ++)• The sets of points on the EDS associated with symbol 
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codes (++; ++) and ( — h; — h) are actually empty. This is due to the fact that, when one 
samples the EDS and tries to solve the equation K(Ix, J 2 ) = E > 0, there is no solution for 
which the two actions are negative simultaneously. 

Figures [I] and |9] show the different parts of the EDS for the two subclasses (/i+; ix~) and 
Different types of trajectories on the EDS are represented by different colors: 
red for type I2, green for type li and blue for type 2. 

The interpretation of these results is quite clear. The fate of trajectories in each of these 
subclasses is apparently controlled to a large extent by the dynamics associated with the 
first saddle plane. Classifying the results according to the symbols in this plane, trajectories 

in the subsubclass (++; H — ) U (++; ) exhibit mainly type I2 behavior, as expected for 

trajectories which either cross the first saddle direction from Q\ < to Qi > or from 
trajectories which remain on the half of the plane (Qi,P\) with Q\ > 0. Referring to Figure 
[7| we see that these trajectories are those for which the bond r 2 is anticipated to break. In 
the same manner, trajectories in the subsubclass ( — h; ) U ( — h; H — ) are mainly type \\. 

An important observation is the fact that type 2 trajectories (those exhibiting double bond 
fission) are not distributed randomly on these parts of the EDS; rather, type 2 trajectories 
actually occur at the boundary between type li and I2 (cf. refs !101[ 1102]) . 

We now turn to the second class of trajectory, consisting of sets (fx— ; U ; h~)- 
Again, this class splits into two subclasses (fx— and ;ii— ), and those subclasses 



split again according to the dynamics in the first saddle plane. Figures 10 and 11 show the 
results for these two subclasses. The situation is now very different from the previous case. 
All the trajectories belonging to these two subclasses enter the well region. The subsequent 
dynamics is then much more complex, and we are unable to predict which bond will be 
broken for a trajectory initiated on this portion of the EDS. 

It should nevertheless be emphasized that the distribution of the different types of trajec- 
tory within these parts of the EDS is by no means random and that we readily recognize a 
pattern of red and green 'stripes', corresponding to alternating single bond breakages. The 
interesting question then arises as to what happens at the boundary between types li and 
I2. In order to answer this question we sampled a line of initial conditions on the EDS 
which intersects these alternation of red and green strips. The results are shown in figure 



12 Figure [I2[a) shows trajectory type versus the P\ coordinate (used to sample the line 



of initial conditions); the results indicate an alternation of the two types li and 1 2 . On 
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the scale of Figure |12|(a), the transition between one type to the other seems to be abrupt; 
however, sampling more densely we get the results shown on Figure pjfb), where now type 
2 behaviour appears at the boundary between type li and 1 2 . 

The interpretation of this result is again very simple. As the energy at which we sample 
the EDS is greater than the energy of the EP4 saddle, trajectories in principle have enough 
energy to pass over the index 2 saddle and break 2 bonds. By continuity, as we pass from 
the set of initial conditions in which the single bond r\ breaks to the set where bond r 2 
breaks, we traverse the region of the EDS where both bonds break (cf. again the work of 
Chesnavich, refs QEH EEJ). 



2. Energy below the index 2 saddle and 'bond healing' 

We now turn to the behavior of trajectories initiated on the EDS at energies below the 
energy of the equilibrium EP4. 

As discussed in the previous subsection, the results are easier to understand if we consider 
separately the part of the EDS consisting of trajectories which immediately escape from the 
well region and the part for which trajectories enter the well. 

As before, the first class of trajectories is i\— ) U(/i+; As for the higher energy 
case, type 2 trajectories appear at the boundary between types li and 1 2 and we obtain 
mainly type li behavior for subsets having symbols (— ; +) and (— ; — ) for motion in saddle 
plane 1 and conversely mainly type 1 2 for subsets having symbols (+; +) and (+; — ). These 



two subclasses are represented on figures [13] and [14 

The second class, as for the case of energies above that of the saddle, show a complex 
dynamics with alternation of type li and 1 2 trajectories. The two subclasses of this class 



are represented on figure 15 and 16 Again, our interest will focus on what happens at the 
boundary. If we sample along a line of initial conditions we see an abrupt change between 
these types; however, even if we sample on a very fine grid we do not see type 2 trajectories 



at the boundary. Figure 17 shows the global representation of the EDS in the toroidal 



representation at energy E = —0.03. 



In Figure 18 we examine the trajectory exit time (the time it takes for a trajectory to 



actually dissociate, either along bond r\ or r 2 ) along the line of initial conditions. Figures 



18 a), (b) show the trajectory type and corresponding exit time, respectively, along the 
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sampling line. At the boundary between types li and I2, the exit time becomes very large, 
and appears to diverge as the sampling density is increased (results not shown here). We 
conclude that, at the boundary between type li and 1 2 trajectories, there exists a set of 
measure zero for which the exit time is infinite. In other words, these are trajectories that are 
trapped in the well region for t — > +00. The interpretation of this result is quite familiar [54J: 
trajectories do not have enough energy to overcome the barrier for double bond breakage so 
that, in order to dissociate, either the bond or the r 2 bond must break. Between these 2 
possibilities, we have trajectories that take an infinite time to 'decide' between dissociation 
channels, and so remain trapped in the well. (Note that, although we initiate trajectories in 
the vicinity of index 2 saddle, trajectories can dissociate either by passing close to the index 
2 saddle or by passing through the 'usual' transition state associated with one of the index 
1 saddles. In the present work, we do not explore the interesting and important dynamical 
interplay between index 1 and index 2 saddles (cf. ref. HOOp . nor do we investigate the 
possibility of defining a 'global' dividing surface encompassing both types of saddle |lUUj .) 

Examination of those trajectories which exhibit a large exit time suggests a connection 
between the form of these trajectories in configuration space and the so-called bond healing 
phenomenon [631 EH EH] ■ 



Figure 19 shows two examples of trajectories having a large exit time. The phenomenon 
of bond healing (as we interpret previous discussions [631 EH EH] of the concept) refers to 
the incipient breakage of one bond of the chain, followed by a recombination or 'healing' 



of the bond rather than dissociation. The trajectories shown in Figure 19 exhibit precisely 
this behavior; the trajectories pass well outside the nominal transition state (index 1 saddle) 
before returning to the well region before ultimately dissociating. In fact, these trajectories 
exhibit what might be called alternating bond healing, where one bond almost breaks and 
then reforms and then the other bond breaks and reforms, and so on and recombined. Such 
trajectories oscillate in an anti-diagonal direction in the (r l5 r 2 ) plane. 

An important task for future investigation is the identification and computation of in- 
variant phase space objects (periodic orbits, NHIMs) which are responsible for trapping 
trajectories in their vicinity, leading to divergent exit times. 
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V. CONCLUSIONS AND PERSPECTIVES 



In this paper we have investigated the fragmentation dynamics of an atomic chain under 
tensile stress. We have analyzed the location and types (indices) of equilibria for the general 
n-particle chain, and have noted the importance of saddle points with index > 1. 

For an n = 2-particle chain under tensile stress the index 2 saddle is of key significance for 
the dynamics. Building upon previous work, we apply normal form theory to analyze phase 
space structure and dynamics in a neighborhood of the index 2 saddle. We define a phase 
dividing surface (DS) that enables us to classify trajectories passing through a neighborhood 
of the saddle point using the values of the integrals associated with the normal form. We 
also generalize our definition of the dividing surface and define an extended dividing surface 
(EDS), which is used to sample and classify trajectories that pass through a phase space 
neighborhood of the index 2 saddle at total energies less than that of the saddle. 

Classical trajectory simulations are used to study single versus double bond breakage for 
the n = 2 chain under tension. Initial conditions for trajectories are obtained by sampling 
the EDS at constant energy. We sample trajectories at fixed energies both above and below 
the energy of the saddle. The fate of trajectories (single versus double bond breakage) is 
explored as a function of the location of the initial condition on the EDS, and connection 
made to the work of Chesnavich on collision-induced dissociation |101[ 1102] . A significant 
finding is that we can readily identify trajectories that exhibit bond healing [63, EH EE]; 
such trajectories pass outside the nominal transition state for single bond dissociation, but 
return to the well region possibly several times before ultimately dissociating. Identification 
of the invariant phase space structures associated with trapped trajectories is a topic for 
future investigation. 
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Appendix A: Sampling the extended dividing surface for index k saddles 



The concept of the extended dividing surface (EDS) (as introduced in Section III) is of 
essential importance for our study. The precise definition of this phase space object is given 
in terms of the normal form Hamiltonian. A procedure for sampling phase space points 
on the dividing surface associated with index-2 saddles was introduced and implemented in 
Ref. 99. In this Appendix we present another method of sampling the Normal Form EDS 
which extends quite naturally to the case of index k saddles where k > 2. 

We start by recalling the definition of the saddle actions in terms of the (Qi, Pi) coordi- 
nates: 

h = \(P?-Ql), 1 = 1,...,*. (Al) 
As the actions I{ are constants of the motion, we have in each saddle plane the equation: 

P 2 O 2 

±*--**± = l, (A2) 

which is just the equation of a doubly-branched hyperbola. There is then a very sim- 
ple parametrization of the hyperbola in terms of one parameter ti, where the form of the 
parametrization depends on the sign of Jj. For the case Zj > we have: 
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In ref |99]the dividing surface was defined as a codimension one surface within the energy 
surface consisting of phase points which minimized the distance from the origin, where the 
distance D was defined in configuration (Qi) space: 

k 

(A5) 



D 



1 k 



i=i 



The use of the distance D in ref. |99| was appropriate because the two parameters R and 
9 used to parametrize the DS specified the location of phase points in the (Qi,Q2) plane 
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(index 2 case). Conservation of energy together with the minimum distance condition then 
served to define the corresponding momenta (P 1; P 2 )- 



If we use the parametrization of the hyperbolas of eq. (A3) and (A4) it is more natural 



to introduce a new distance defined in the full phase space (cf. eq. (3.12)): 



i=i 



(A6) 



The derivatives with respect to time of the two distances D and D are in fact identical, so 
that for both cases the minimum distance condition is: 

k 



(A7) 



Independent of the sign of Ii, the following equation then defines the constant energy EDS 
(K(I l7 ...,I n ) = E): 

(A8) 



i=i 



cos(tj 



In terms of the variables Xi = tan(tj)/cos(tj), we have: 



Y,2Ai\Ii\X t = 0. 

i=i 



(A9) 



As noted above, the EDS is composed of 4 fc parts for the case of an index k saddle. For 
each saddle plane, the hyperbola has four branches, two with /, > and two with Jj < 0. The 
four branches are in 1:1 correspondence with the symbolic notation introduced in Section 
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From the equation Xi = tan(tj)/cos(tj) we obtain: 
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which yields 

z] - Zl -X? = (A12) 



in terms of the new variable Zj = l/cos 2 (tj). The discriminant A = 1 + AXf of eq. (A12) is 
always positive and there are two roots: 



/ x 1±VA . . 

= — ^— ■ (A13) 

However, z is also a square (z = cos(tj)" 2 ) and the solution is always negative or zero. 
So, the only possible solution is actually (zj)+, so we finally have for ti the expression: 



U = cos 1 



±1 



i + Va 



(A14) 



If we set a = y i+v^ ' we see < a < 1 and — 1 < — a < and consequently 
cos -1 (a) G [0; |] and cos -1 (—a) G [| ;7r]. It is straightforward to choose the appropriate root 
depending on which branch of the hyperbola we are working on and making the appropriate 
translation of the set [0; ir) to either [— |; |] or [|; 

To sample the EDS, we first have to fix the value of the energy E, which can be either 



positive or negative. We must therefore find a set of actions Jj for which eq. (3.11) holds. 
This can be done by choosing n — 1 actions and solving numerically for the value of the last 
action I s . This set of actions can be chosen in a systematic way. We can specify maximal 
or minimal values for each action depending on the part of the EDS we are sampling and 
sample systematically from zero to the maximal or the minimal value of each actions. The 
determination of these maximal or minimal values of the actions is related to the accuracy 



of the normal form Hamiltonian and its associated coordinates transforms (3.2). The check 
of the accuracy of the normal form Hamiltonian is an important topic and we discuss some 
aspects of this topic relevant to this work in the next paragraph. 

Normal form theory is procedure which enable us to construct a "simple" Hamiltonian 
which approximates the dynamics of a given Hamiltonian in a neighbourhood of an equilib- 
rium point of this Hamiltonian. The outputs of this procedure are a normal form Hamilto- 
nian and two sets of coordinates transformations which enable us to navigate between the 
two sets of physical and normal form coordinates. The construction of the normal form 
Hamiltonian relies on a Taylor expansion about an equilibrium point and a normalisation 
process which results in the normal form being valid only in a certain neighbourhood of 
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the equilibrium point under consideration. Therefore we expect the normal form to "break 
down", i.e. not be an accurate approximation of the physical Hamiltonian, when we leave 
this neighbourhood. The major question here is how to define this neighbourhood or, in 
other words, how to measure the accuracy of the normal form? For Hamiltonian systems 
there is a quite natural measure of the accuracy of the normal form which is the conserva- 
tion of energy. The idea used here is to compare the energy provided by the normal form 
Hamiltonian with the energy provided by the initial Hamiltonian and decide a threshold for 
energy "disagreement" between the two: 



\E mi -E NF \ < e. (A15) 
For phase space points in normal form coordinates we compute the energy with the normal 



form Hamiltonian and use the backward transformations (3.2) provided by the normalisation 
process to obtain phase space points in physical coordinates with which we can compute 
the energy using the physical Hamiltonian. Starting from phase space points "close" to 
the equilibrium point, and increasing the distance from this equilibrium point gradually, we 



can construct a connected set of phase space points for which (A15) holds, and which for 
the purpose of this work defines the neighbourhood of validity of the normal form. This 
neighbourhood, denoted C e , depends on e. The parameter e actually quantifies the "error 
threshold" in the normal form. For this work the value of e was taken as a factor 10~ 3 times 
the physical energy of the energy surface on which the EDS was sampled. 

We now return to the problem of computing the maximal or minimal values of the actions. 
We can differentiate between the case of a saddle mode action and a bath mode action. For 
a bath mode with action /& the dynamics is confined to a circle of radius r = y/2Ib in this 
bath mode plane. To determine the maximal value of h we set all other coordinates related 
to the other modes (bath and saddle modes) to zero so that the phase space point is located 
at the origin of all the mode planes except the one we try to determine the maximum value 
of the action. For a certain value of the action we calculate the energy e-agreement for 
samples all around the circle of radius r = As we increase the action If, the radius 

of the circle will increase and we will find some points on the circle which do not belong to 
C e . So the maximal value of the action lb will be the one corresponding to the last circle 
for which all the points on the circle belong to C e . For the action corresponding to a saddle 
mode we can have two cases depending on which part of the EDS we are sampling but the 
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procedure is quite identical for the two cases. For a saddle mode of fixed action I s on a 
certain part of the EDS the dynamics is constrained on one branch of a hyperbola. As for 
the case of a bath mode we set all coordinates related to the other modes to zero. On the 
branch of the hyperbola the minimal distance from the origin in the saddle mode plane is at 
t = or 7r depending on the part of the EDS. So if we fix t to or it and increase (I s > 0) 
or decrease (I s < 0) we will pass a point for which the energy e-agreement between normal 
form Hamiltonian energy and initial energy does not holds which means that we have left 
C t . This determines the maximal or minimal value of the action I s . 

Having selected a particular set of bath mode and saddle actions we sample points which 
belong to the EDS by assigning values to the the variables that are complementary to the 
actions. For each of the bath mode planes we sample the angle variable conjugate to the 
actions. If a bath mode action is lb, a conserved quantity, we have h = \{ql +pl)- The 
angle 8, < 9 < 2tt, parameterizes a circle of radius r = y/2If, in the {qb,Pb) plane. 

For the saddle planes, we need to solve the minimum distance equation to sample points 



on the EDS. Equation (A9) provides us with a very simple relation among the variables X. 



i ) 



i = 1, . . . , k. By fixing all of those variables except X\ and X2, for example, equation (A9) 
describes the equation of a straight line with slope — xfjrfj an d F-intercept £)* =3 ~ xnjn Xj. 
For the case of only two saddle planes, k — 2, it reduces to the equation of a straight line 
passing through the origin. After determining the Xi, we compute the ti as explained above 
and finally the (Qi, Pi). 
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FIG. 1: (a) Potential function of eq. (2.1 ) for different values of the force parameter /. (b) Locations 



of the equilibrium points of the 1 dof Hamiltonian (2.1) as a function of the force parameter /. 
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FIG. 2: Definition of 'external' (space-fixed) coordinates (xi,X2) and 'internal' (bond) coordinates 
( r i, r 2) for the n = 2 atom Morse chain. 
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FIG. 3: (a) Contour plot of the potential energy surface for the n = 2 Morse chain with / 
(b) Energies of the different equilibria as a function of the force parameter /. 
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FIG. 4: Schematic representation of the dynamics in the two saddle planes in physical normal form 
coordinates. 
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FIG. 6: Toroidal representation of the extended dividing surface. 
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FIG. 7: Potential energy contours and projections of eigenvectors obtained by linearizing Hamil- 
ton's equations about the equilibrium point EP4. 
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FIG. 9: Subsets of the EDS (/i+;»i+). (a) (++;-+). (b) (-+;++). 
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FIG. 10: Subsets of EDS (fi-,k+). (a) (+-;++). (b) (+-;-+). (c) ( ;++). (d) ( ;-+). 
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FIG. 12: (a) Trajectory type as a function of location along the ID cut of the EDS: Q\ = 0.2, 
Pi = —0.35-0.0, E = 0.03. (b) Magnified segment showing the appearance of type 2 trajectories 
at the boundary between type li and I2. 
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FIG. 13: EDS subsets (a) (++;+-). (b) (-+; ). 



48 



2 




' (a) 




' (b) 


n 
4 










t 2 o 

n 




"\ . 


■■ i - ' if mil—" - 




"4 












— h — i 


1 

(c) 


1 


1 — t — 


n 
4 


- 1 








?2 








1 


n 
"4 

;r 











- - — £ 3^ 5g ~in 

"2"4 U 424 42 



FIG. 14: EDS subsets (a) (++;-+). (b) (-+;++). (c) (++;++). (d) (-+ 



49 




FIG. 15: EDS subsets (a) (+-;++). (b) ( ; 
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FIG. 17: Global representation of the EDS using the toroidal representation. Energy E = —0.03. 
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FIG. 18: (a) Trajectory type as a function of location along the ID cut of the EDS: Q\ = 0, Pi 
—0.1-0.1, E = — 0.03. (b) Trajectory exit time s in time units along the cut. 
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FIG. 19: Examples of trajectories exhibiting the phenomenon of 'bond healing'. 



54 



